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New algorithm for quantum state estimation based on the 
maximum likelihood estimation is proposed. Existing tech- 
niques for state reconstruction based on the inversion of mea- 
sured data are shown to be overestimated since they do not 
guarantee the positive definiteness of the reconstructed den- 
sity matrix. 



State reconstruction belongs to the topical problems of 
contemporary quantum theory. This sophisticated tech- 
nique is trying to determine the maximum amount of 
information about the system-its quantum state. Even 
if the history of the problem may be traced back to the 
early days of quantum mechanics, till quantum optics 
opened the new era of the state reconstruction. The- 
oretical prediction of Vogel and Risken was closely 
followed by the experimental realization of the suggested 
algorithm by Smithey et. al. j^]. Since that time many 
improvements and new techniques have been proposed 
|p|-p2[, to cite without requirements to completeness at 
least some titles from the existing literature . Even if 
the method comes from optics, similar methods are cur- 
rently being used also in atomic physics as quantum en- 
doscopy [Q. Homodyne detection of quadrature opera- 
tor with varying phase of local oscillator {x^, (j>) was used 
as the measurement in the original proposal The 
algorithm served for determination of the Wigner func- 
tion W{x,p) and also other quasiprobabilities represent- 
ing the density matrix. Measurement of rotated quadra- 
ture operator may also be used for direct evaluation of the 
coefficient of density matrix in number state representa- 
tion pm,n and for the analysis of multimode fields [Q. 
Simultaneous measurement of the pair of quadrature op- 
erators {x,p) using double homodyne or heterodyne de- 
tection yields directly the Q-function Q{a) ||]. Surpris- 
ingly easy technique was suggested by Walentovitz and 
Vogel and by Banaszek and Wodkiewicz |0 . Mixing 
of the signal and coherent fields with controlled ampli- 
tude on the beam splitter may serve for reconstruction 
of the Wigner function and other distributon functions 
using the photon counting only. Techniques similar to 
the quantum state reconstruction have been suggested 
for indirect observations of particle number, see for ex- 
ample Ref. ||l5|. Though the techniques are different as 
far as the practical realization is concerned, they all may 
be comfortably represented by the formalism of general- 
ized measurement [|6| . As well known, any measurement 
may be described using the probability operator measure 



(POM) n(^), being any positively defined resolution of 
identity operator n(^) > 0, J d^U{^) = i. Probability 
distribution of the outcome predicted by quantum the- 
ory is 



,(e)=Tr[/}n(a], 



(1) 



where p is the density matrix of the state. The measured 
variable ^ represents formally the registered data being in 
general a multidimensional vector with the components 
belonging to both the discrete and continuous spectrum 
as shown in above mentioned examples. The key point 
of the existing reconstruction techniques-inversion of the 
relation (|l|) , represents nontrivial problem. Solution may 
be formally written as an analytical identity 
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Wp{a) being a representation of density matrix. In or- 
der to find the representation M^(a) of a density matrix 
corresponding to an unknown signal, the existing recon- 
struction techniques apply the relation (||) on the actually 
detected statistics w{^). 

Apart from the fact how ingeniously the individual in- 
versions have been done, this treatment is essentially im- 
proper for the application in quantum theory. Partic- 
ularly, it may represent a density matrix only for such 
measurable probability distributions, which are given ex- 
actly by the relation (|^) . Deviations between actually de- 
tected w{£,) and the true statistics Wp{^) are not allowed, 
since they may spoil the positivity of reconstructed den- 
sity matrix. This algorithm anticipates therefore the ab- 
solute precision impossible in quantum theory. There are 
at least the following imperfections of detected statistics 
w{^), which should be taken into account: i) Sampling 
error caused by the limited number of available scanned 
positions of continuous variable at which the measure- 
ment was done, ii) Counting error caused by limited set 
of available data counted at each position. For exam- 
ple, in Ref. Q the former one is caused by the division 
of quadrature x^ into 64 bins and phase into 27 values, 
whereas the later one by the detection of quadrature x^ 
at each bin. Other errors such as imperfections of detec- 
tors or external noises, may appear in practise as well. In 
quantum case the algorithm based on inversion provides 
a result, but does not guarantee the positive definiteness 
of the reconstructed density matrix 1 1^ . In the example 
of Ref. j|] the positive definiteness of reconstructed ma- 
trix has not been checked explicitly, but can be judged 
according to the papers flp|Jl8|. Here the negative part 
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of photocount distribution indicates the spoihng of pos- 
itive definiteness. Even if there is a connection between 
the dimension of Hilbert space where this happens and 
the number of phases |l8| , a righorous way how to treat 
the positive definiteness within reconstruction has not yet 
been suggested. As pointed out by Jones jl9| in his Ref. 
12, the failure of similar methods is the rule rather than 
the exception in the case where the data underdetermine 
the solution. This happens also in the case of more dense 
data, even if the region of ill behaviour is shifted to less 
obvious manifestations. These methods are considered 
as satisfactory only because some additional information 
in the form of data smoothing is used and mathemati- 
cal difficulties are neglected. Instead of inversion of the 
detected data, a technique motivated by quantum infor- 
mation theory (l^j2^ and by phase shift estimation , 
will be suggested in this Rapid Communication. Previ- 
ous reconstruction techniques will be embedded into the 
common scheme based on maximum likelihood estima- 
tion. 

Many parameters characterizing the quantum state 
should be estimated in state reconstruction. As pointed 
out by Helstrom [l^ this may be done restricting the 
dimension of Hilbert space, and accepting some resid- 
ual uncertainty. Similarly, Jones investigated the 
fundamental limitations of quantum state measurement 
using Bayesian methodology. On the contrary, the real- 
istic measurements as in the existing techniques will be 
anticipated here. Assuming the repeated (or multiple) 
measurement performed on the n copies of the system, 
the output of observation may be parametrized by the set 
of states (projectors) formally denoted as |yi), . . . , |j/„), 
repetition of a particular outcome being allowed. Pure 
states represent here the case of sharp measurement, 
whereas unsharp measurement involving the finite res- 
olution should be represented by an appropriate POM. 
Since formal considerations are valid for both these cases, 
the notation of sharp measurement will be kept in the 
following for the sake of simplicity. Maximum likelihood 
estimation ascribes to such a measurement the state p 
maximizing the likelihood functional 
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The aim of this contribution is to find this state and to 
clarify the fiuctuations of such a prediction. As the math- 
ematical tool, the inequality between the geometric and 
arithmetic averages of non-negative numbers qi will be 

used (n"?*)^^" — nY^l 'ii- The equality is achieved if 
and only if all the numbers qi are equal. The variables 
will be formally replaced by qi = Xi/ai, where Xi > Q 
are positive and > are auxiliary positive nonzero 
numbers. In the following the n dimensial vectors will be 
denoted by boldface as a, x, y, etc.. Assume now that 
the numbers qi are chosen from the given set of values 
so that the value qi appears ki times in the collection 
of n data. Hence ki represents the frequency, fi — ki/n 



being the relative frequency fi = 1. Parametrization 
revealing explicitly the frequency will be denoted by up- 
per prime in sums and products, indicating that index 
runs over spectrum of different values. Without loss of 
generality the variable x may be interpreted as proba- 
bility ^[ Xi — 1, since the normalization may always be 
involved in auxiliary variables a. The relation, known as 
Jensen's inequality ||22| , then reads 



n; 



< 



f ^ 1. 
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In this form it represents remarkably powerful! relation 
since the equality sign may be achieved for an arbi- 
traryprobability x = a. For example, the Gibbs inequal- 
ity |20|] follows as a special case chosing the parameters 
— fii since the inequality (|^) may be rewritten as 
~ Yli /i 111 ~ ^ 0- These formal manipulations are tightly 
connected to the maximization of likelihood function. 
Using the definition 
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tti being a subject of further considerations, the likelihood 
functional may be simply estimated as 



mf-^]X\{yM 



< 



n.af Tr{pi?(y,a)}. (6) 



The operator R is given, in general by nonorthogonal, 
decomposition as 



i?(y,a)=^'^|y,)(y, 



(7) 



Relation (^ simply follows from the definition (^ and 
from the inequality (^. Further treatment is distin- 
guished by the following specification of auxiliary param- 
eters a: 

Reconstructions of wave function 
Condition Oi — fi tends to considerable simplifications. 
Since the measurement need not be complete R{y, a = 
f) < i, the right-hand side of the relation (ph reads 



(6)=n'//^Tr{pr.|y.)(y.|}<n'/: 



(8) 



This represents a state-independent upper bound. The 
necessary condition for equality sign in (|^) is given by 
the conditions {yi\p\yi) /oi — const for any i, whereas 
the equality sign appears in relation (|^) for complete 
measurements. These relations together with the nor- 
malization of relative frequencies tend to the necessary 
condition for searched state p 
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This is nothing else as the experimental counterpart of 
the relation (|l]) and hence the starting point of recon- 
struction based on inversion. The relation (|) may be 
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simply inverted in the case of orthogonal measurements, 
which may be considerd as complete on the given sub- 
space, tending to the solution 



Pf 



(10) 



Unfortunately, such measurements do not reveal infor- 
mation about full density matrix since the nondiagonal 
elements are lost, as for example in the case of particle 
number measurement. Techniques dealing with orthogo- 
nal measurements are therefore not suitable for full state 
reconstruction, which should be based on the usage of 
nonorthogonal states. On the other hand, in these cases 
the completeness and the existence of a solution of the 
equation (^ cannot be guaranteed. Quantum analogy 
of Gibbs inequality corresponds to overestimated upper 
bound and tends to the conditions imposed by recon- 
struction techniques. 

Maximum likelihood estimation 
The problems with existence of a state achieving the up- 
per bound descends obviously from the fixing of the aux- 
iliary parameters a. The remedy is to keep them free as 
a subject of further optimization. For any positively de- 
fined operator B — '^j^Wbi){bi\ and density operator p 
the simple Lemma holds 



Tt(pB) < max A; 



(11) 



The quality sign is achieved for density matrix corre- 
sponding to the spectral projector of operator B with 
maximal eigenvalue. Using this Lemma the estimation 
of the right hand side of the inequality (0) than reads 



(6) <A(y,a)[]' 



Ji 
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where A(y, a) denotes formally the maximal eigenvalue of 
the operator R{y, a) with the corresponding eigenvector 
|'0(y,a)). Equality signs in the chain of inequalities are 
achieved simultaneously if and only if 



— const. 
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independently on the index i. Finally, maximum likeli- 
hood estimation determines the desired state as \ip{y, a)), 
where vector a is given by the solution of the set of non- 
linear equations ([l^) . The uncertainty of such a quantum 
state estimation may be, according to the Bayesian for- 
mulation 1^ , characterized by the likelihood functional 
(^. Since the interpretation of the probability distri- 
bution on the space of states is rather complicated, the 
uncertainty of the prediction may be involved in an alter- 
native way. The measured data are fiuctuating according 
to the distribution function P(y) depending on the true 
state of the system. Fluctuations of quantum state es- 
timates may be represented by the sum of independent 
contributions 
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This density matrix shows how closely the maximum like- 
lihood method allows to estimate an unknown state hid- 
den in the measured statistics P{y)- Unfortunately, the 
proposed method is rather complicated and examples of 
reconstructions specified above should be solved sepa- 
rately case by case. Considerable technical difficulties 
may be caused, for example, by possible degeneracy of 
operator B refiecting the structure of performed quan- 
tum measurement. This particular questions are beyond 
the scope of this contribution and represent an advanced 
program for further re-interpretation of existing recon- 
struction techniques. 

Developed technique may be illustrated on simple but 
theoretically worth examples. Quantum state recon- 
struction after the measurement of a Hermitian operator 
with orthogonal spectrum is the simplest problem. Solu- 
tion corresponds to the application of Gibbs inequality, 
since the relation (^ may be solved in this case. Quan- 
tum state is then reconstructed after each measurement 
by the density matrix (^0|) . This is a consequence of the 
possible degeneracy of the operator B mentioned above. 
The treatment based on the Gibbs inequality is overesti- 
mated in general. Provided that (^ is fulfilled in some 
special cases, then the solution should coincide with the 
prediction of maximum likelihood estimation. 

Simple are also the cases of strongly underdetermined 
data, when the state is estimated after single detection 
n — 1. Assume for concretness the standard "measure- 
ment of Q-function" corresponding to detection of co- 
herent states \y) — e^°^^^ '^jO). If the value y is detected, 
the system is with the highest likelihood just in the state 
\y). Provided that system was in coherent state \a), the 
output fluctuates as \{a\y)\'^ /tt. Estimation after single 
detection then yields the density matrix of superposition 
of coherent signal a and the thermal noise [ p3[ with mean 
number of particles equal to 1, 
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Difference between the true state and its estimation is 
negligible in the case of classical fields, but considerable 
in quantum domain. 

Estimating the quantum state after multiple detection 
of coherent states, the matrix R should be diagonalized. 
Using the assumption for eigenstates as If) — Vi\yi), 
linear equations for desired coefficients Vi and eigenvalues 
A follow as 



fk 
CLk 



(15) 



where Cui = C*^, = {yk\Ui),Cii = 1. This solution de- 
termines the coefficients a according to the relation ( |T^ ) 
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^ lYliViCki]'^ /ak = const for any index k. Let us il- 
lustrate this strategy on the case of double detection 
n = 2 yielding the values t/i and t/2- Parameters are 
given as /i = /2 = 1/2 and without loss of general- 
ity oi = l,ai/a2 — X. The secular equation for maxi- 
mal eigenvalue A reads — (1 + a;)A + x — a;|Ci2p = 0, 
yielding easily solutions for maximal eigenvalue and its 
eigenvector. Equations ( |l3|) impose single condition as 
|Ci2|A = y/x{X — 1 -I- |Ci7f ). This nonhnear system of 
equations may be easily solved yielding expected solu- 
tion as A = 1 + |Ci2|,x — 1. Projector is given by the 
normalized Schrodinger cat-like state 

Density matrix "reconstructing" the coherent state is 
then given as 

PMLE = ^ y"d2yid2y2e-|^^-"l'-|^^-"lV>('^|- 

Proposed method describes easily the cases, where the 
data seem to be underdetermined. There is also a strong 
effort to apply the developed technique to the case of 
large data sets estimating properly the quantum state in 
the cases of realistic measurements. 

Even if the problem of positive definiteness used for 
motivation may seem as nit-picking, it has far-reaching 
consequeces. Reconstruction methods based on inversion 
of detected data are valid only if complete information is 
available and fail, if the information is limited by quan- 
tum theory. The method based on maximum likelihood 
suggests, how to treat the state reconstruction in this 
quantum domain. Since quantum state comprises max- 
imum possible information about the system, its proper 
description is of fundamental interest. For example, the 
strategy of indirect measurement observing primarily the 
wave function and deriving all information about desired 
variable consequently, represents a general scheme for an 
universal "measurement of everything" . The fee paid for 
such an observation should be obviously the accuracy of 
the detection of a particular observable, since the obser- 
vation of wave function involves the registration of non- 
commuting variable, too. The fundamental distinction 
between universal and accurate measurements may eas- 
ily disappear provided that an improper description of 
quantum objects is used. Particularly, this is just the 
case of existing reconstruction techniques, where "nega- 
tive probabilities" may appear as a consequence of inad- 
equate semiclassical treatment. 
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